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Abstract 


Vi  X 

{/ 

aaveslan  nonparametr ic  estimators  ot  the  survival  function,  the  failure 
rale  function,  and  the  densitv  function  are  obtained  using  Jump  processes  as 
prior  distributions  on  the  space  of  increasing  failure  rate  functions.  Th« 

’ ur  processes  are  more  appealing  intuitively  than  previously  used  Dirichlet 
proi  esses,  integral  gamn.i  process,  and  processes  neuttal  to  the  right  and  hav«  a 
meaningful  phvsital  interpretation.  Examples  are  given  and  the  estimates  are 
compared  with  tin  maximum  likelihood  estimates.  In  addition,  the  Bavesian 
nonparamet r ie  estimators  are  presented  for  arbitrarily  r 1 ght -censor ed  observations 


f  Accc  r 


Kc^-  Words:  Bayesian  nonparamet r ic  estimation;  Jump  process  priors;  Life  testing; 
Failure  rate  function;  Right  censorship;  Reliability. 


1.  INTRODUCTION 


The  problem,  of  Bayesian  nonpararae trie  estimation  of  a  survival  or  rellabil 
ltv  function  has  been  considered  in  recent  years  bv  several  authors  (for 
example,  see  Ferguson  '.1973,  197*.),  Ferguson  and  Phadia  (1979),  Lo  (1979), 
and  Susarla  and  Van  Rw  in  (1976)).  These  authors  have  utilized  malulv  the 
Dirlchlet  processes,  integral  gamma  pr  'cesses,  or  processes  neutral  to  the 
right  a>-  prior  probability  distributions  over  the  space  of  all  distribution 
functions  of  interest.  The  Dirlchlet  process  priors  put  all  of  the  probability 
on  discrete  distributions  while  process*-  neutral  to  the  right  and  integral 
gammu  processes  are  somewhat  difficult  to  use  in  practice.  Also,  if  further 
prior  inf  orr.a  t  Ion  concerning  the  distribution  function  F  to  be  estimated  is 
known,  the  Dirlchlet  process  priors  mav  not  be  entirely  appropriate.  For 
example.  If  It  Is  assumed  that  F  has  an  increasing  failure  rate  (IFF) 
function,  the  Dirlchlet  process  would  not  be  suitable  since  F  is  absolutely 
continuous  on  at  least  part  of  its  support.  In  this  paper  we  propose  a  prior 
process  on  the  failure  rate  function  for  Bayes  estimation  of  the  failure  rat< 
and  survival  functions  under  the  assumption  of  IFK.  The  proposed  prior  is 
very  practical  and  has  a  physical  interpretation  similar  to  that  of  a  shock 
model  with  shocks  occurring  as  a  Poisson  process.  No  assumptions  are  made 
about  F  except  that  it  is  IFF. 

Let  r(t  denote  the  true  failure  rate  function  of  the  unknown  distribu¬ 
tion  function  F  and  let  F(t)  denote  survival  probability  at  t  2  0  .  We 
assume  that  F  has  a  density  function  f  <  t )  “  r(t)F(t)  .  Nonparamctr ic 
estimation  of  the  failure  rate  function  has  been  considered  by  Watson  and 
Leadbetter  (1964a,  1964b)  and  bv  Barlow  and  van  Zwet  (1971),  among  others. 

In  particular,  Grenander  (1956)  and  Marshall  and  Proschan  (1965)  obtained  the 


maximum  like.ihood  estimator  of  r(t)  under  increasing  assumptions.  Their 


estimator  was  a  step  function  and  took  the  value  ♦»  at  the  largest  value  of 
the  sampl«  .  f'adgetl  anu  Wei  (1979)  have  recently  obtained  the  maximum  like¬ 
lihood  estimator  of  monotone  failure  rate  functions  for  arbitrarily  right- 
censored  survival  data. 

The  prior  probability  distribution  on  r(t)  which  we  propose  is  moti¬ 
vated  by  the  fact  that  r(t)  can  be  closely  approximated  and  estimated  by  a 
step  function.  We  assume  throughout  the  paper  that  r(t)  is  increasing 
(meaning  nondecreasing)  for  t  >  0  so  F  is  an  1FK  life  distribution. 

Suppose  the  experimenter  has  prior  information  that  he  can  quantify  as  a  jump 
process  for  rtt)  with  constant  jumps  of  si2e  t  at  times  T  >  0  ,  i  • 

1,2,  •••  ,  where  T ^  arc  arrival  times  for  a  Poisson  process  {N(t)  :tsO' 

with  intensity  rate  v  .  This  constant  jump  process  thus  describes  the  prior 

probability  distribution  over  the  parameter  space  G  -  (r  :  r  is  an 

increasing  failure  rate  function)  .  This  prior  distribution  over  0  also 

has  the  following  intuitively  appealing  physical  interpretation:  Suppose 

shocks  occur  randomly  in  t iae  to  a  system  or  component  causing  a  certain 

amount  of  damage  to  it  which  increases  the  failure  rate  of  such  components  by 

a  certain  constant  positive  amount  •  .  Assume  that  the  damage  is  cumulative 

and  that  the  shocks  occur  according  to  a  Poisson  process  { N ( t )  :  t  2 0)  with 

intensity  v  .  Then  the  failure  rate  function  can  be  considered,  at  least 

a  priori .  as  a  constant  Jump  process.  This  physical  interpretation  is  similar 

t->  those  probability  models  proposed  bv  Gaver  (1963),  Antelm.m  and  Savage  (196b).  and 

Reynolds  and  Savage  (1^71)  in  which  the  failure  rate  function  was  stochastic. 

In  Section  2  we  will  utilise  the  constant  jump  process  as  a  prior  distri¬ 
bution  on  increasing  failure  rate  functions  to  obtain  a  Bayes  nonparametr le 


estimator  of  r(t)  and  of  the  reliability  or  survival  function  F(t)  ■ 

exp {-/ 1  rtu)du>  based  on  a  random  sampl.  of  size  n  from  F(t)  .  This 
0 

estimator  of  rlt)  turns  out  to  ‘  e  a  smooth  curve  rather  than  a  step  function. 

In  Section  3,  an  example  is  given  using  data  simulated  by  Dubey  (1967)  from  a 
Wcibull  distribution  with  increasing  failure  rate.  Finally,  in  Section  4,  it 
is  shown  that  the  same  techniques  may  be  eisilv  applied  to  obtain  bayeslan 
nonparametr ic  estimators  of  rtt)  and  F(t)  for  arbitrarily  right-censored 
data.  This  problem  was  studied  bv  Susarl.i  nd  Van  kvzin  (1976)  for  Dirichlet 
process  priors  and  by  Ferguson  and  t’hadia  (1979)  for  processes  neutral  to  the 
right . 

2.  THi:  BAYESIAN  NONPAKAMET  It  1ST  I MATO PS  FOR  COMPLETE  SAMPLES 

Let  R(t)  denote  the  failure  rate  function  with  prior  distribution 
given  by  the  constant  Jump  process  described  in  Section  1.  Then  from  the  Jump 
size  <  >  0  and  Poisson  process  {N(t)  :  t  0  .  we  have  R(t)  "  <  N(t)  with 

the  corresponding  survival  function  given  bv 

t  N(t ) 

F(t)  ■  exp.-/  R(u)du  *  exp  ->  (t-T  )  ,  t  i  0  ,  (2.1) 

0  i-1  1 

where  T^  ,  1  *  1,2,  •••  ,N(t)  ,  denote  the  "arrival"  times  of  Jumps  of  size 

*  in  (O.t  !  .  Given  N(t)  ”  k  and  t  ,  It  is  well-known  that  Tj ,  •••  ,T^ 

are  independent  uniformly  distributed  i andom  variables  on  (O.t 

Let  x,  ,  •••  ,x  be  a  given  sample  of  size  n  from  the  unknown  distribu- 
1  n 

tlon  function  F  .  It  is  desired  to  obtain  Baves  estimators  F(t)  and  R(t) 
of  F(t)  and  r(t)  ,  respectively,  with  a  squared-error  loss  function  of  the 
form  L(F,F)  -  / ^  t  F(t)-F(t)  "dw(t)  ,  where  w(t)  is  an  arbitrary  non¬ 


random  weight  function. 


The  likelihood  function  of  the  sample  is 


‘-NUj  .•••.NUn)‘ET|NUJ) i-l,---,n  ^  ki  eXp('*  jf/VV*''  • 


(2.4  > 


where  J  "  (Tj.  •••  .T^  ^  .  The  exponent  in  (2.4)  may  be  rewritten  in  the 

n 

forr. 


k 

n  i 


<  I  l  <*t-T  )  -  -«lk  l  x  +  (k  -k  )  l  x.  +••*+  (k  -k  )* 

1-1  j-1  1  J  1  i-1  1  2  1  l-l  i  n  n-1  n 


] 


1 

♦  <  i  n  l  T  +  (n- 

J-l  J 


i)  I  ♦•••♦  I  V  • 

J’V1  ^n-l*1  J 


Sow,  conditioned  on  k  jumps  in  the  interval  (0,x  J  ,  T  .,  •••  ,T 

•  —  1  *  i  x  ^  |  i  I*  ( 

are  order  statistics  fro®  independent  uniformly  distributed  rando®  variables 
on  (x  ,x  3  •  i  -  2,  **•  ,n  and  T  ,  •••  ,T  are  order  statistics  fro® 


i-1’  i 


ll 


independent  uniformly  distributed  rando®  variables  on  (0,Xj].  Hence,  letting 

Xp  -  0  and  k(1  •  0  .  we  have  for  the  1th  interval  (x^  j.x^]  ,  i  -  1,  •••  ,n 

k  k 

xi  xi  1  Ki 

j  •••  /  exp[  (n-i+1 )  r  "  t  IT  dt  /(x.-x.  .) 

*->  ‘-i  j-.,.,..  i  j-k,.,*.  J 


-  UexpU  (n-i  +  Dx^  -  exp(f  (n-i+1  )Xj  j  )  ]/t  C  (n-i+1)  (xj-x^j  )  ]) 


k  -k 

i  i-1 


By  the  mean  value  theorem,  there  is  an  *1  c(Xj  x  )  for  which 


expi (n-1+1 )  f  Xj  -  exp[ (n-i+1)  cx^  ^ 


(n-i+1 )<x  -x  ) 


-  €  exp[ (n-1+1 )  f  , 


so  (2.5)  can  be  simplified  to  exp'  «  (n-i  +  1  )x^  ( k  ^  — k  ^  j).)  .  Thus,  summing  over 


the  appropriate  Poisson  variables,  the  expectation  in  the  right-hand  side  cf 
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.-)  becomes 


n  -v*r 
t  e 


k  -0  k  -0 
n  n-1 


2  n  n  . 

••  l  (  H  k  )  exp{  <  l  [(n-j+l)x 

kj-0  1-1  j-1 


k  n 


I.Kk-k  ))vn  »(«-«  ,) 
i-j  J  J  1  1-1  1  1  1 


Wi 


/  II  (k  -k  )' 
'  1  i  i-r 

i-i 


(2.7) 


Equation  (2.6)  can  be  solved  for 


■  t  -  1  J»’.  exp  (  (n-i+1 )  <  x  j )  -  exp  ( (n-l*l)(  Kj  j )  ] 


-  ie  (n-i>l)<  (x^-x^j  )  ])/  c  (n-i+1) 


which  can  be  substituted  Into  (2.7)  to  obtain  the  exact  expression  for  the 
denominator . 

To  calculate  the  numerator  of  (2.3)  the  technique  is  the  same  as  for  the 
denominator  except  that  we  have  an  extra  term.  The  numerator  is 

n 

/  F(t)  f(x  >dP 
i-1 

k 

t 

ES<x,>.  •••  .N(xn).S(t)  ET|N(x1)-k1.  •••  ,N(t)-kt  CXp  ‘<  j£j(t'V 


n  “i 


«n  n  l  exp  -  .  [  [  (x.-T  ) 

i-1  1-1  j-1  1  J 


(2.8) 


Note  tnat  the  order  of  summation  corresponding  to  the  first  expectation  in 

(2.8)  depends  on  the  value  of  t  .  If  x(  ^  <  t  $  ,  i-1,  •••  ,n *  1  , 

where  x  -  0  and  x  -  -  ,  then  the  value  of  N(t)  is  between  S(x,  , ) 
u  n* 1  i-j 

and  Nfx^)  wltli  probability  one,  with  N(0)  -  0  and  N(->)  -  »•  .  Hence, 


we  can  simply  relabel  the  x.  s  and  t  as  n+1  new  points  by  defining 


.  J  S  i  -  1 

Vj  ■  >j(t)  -  .  t  ,  J  -  1  .  i  •  1,2,  •••  ,n  . 


xj-l  *  J  '  1  +  1 


Then  similar  to  the  denominator,  the  right-hand  side  of  (2.8)  becomes 


•Wp,  •••  .K<yiiil)iEI:s(y1)-k1 . lj 


n+1 

*  exp!  - c  /  l  (y  -T  ) ]  ■ ) 
j-1  m-1  J 


n*l  n*l  j-1 

iri 


V  k 

n*l  2  n+1  k 


•  t  exp(-vy 


,  i  T  l  ...  I  (  n  k  >v 

W°  k„-°  li*°  J" 


n  1  krki-i 

n  (y  -y  .)  J  J  ex; 

j-1  J  J 


p[  <  l  ((n-j+2)y°  -  l  y  )(k  -k  ' 
j-1  J  m-j  J  J  1 


/  n  (k  -k  ) !  , 
J-1  1  J 


Xl-1  *  1  5  X1  *  1  ’  l'  *"  *n+  1  • 


Therefore,  the  desired  Raves  estimator  F(t)  of  survival  probability 


F(t>  is  given  bv  (2.10)  divided  by  (2.7).  This  estimator  is  readily  comput 


able  bv  electronic  computer  as  the  example  in  Section  3  shows. 


8 


Now.  for  the  bavesian  nonparasietric  estimator  of  r(t)  ,  the  failure 

rate  tuncUon,  wi-  simply  calculate  the  required  numerator  expression  for 

the  posterior  mean  of  K(t)  •  <N(l)  ,  which  is 

n  n  x . 

/  K  ( t )  Hflx^dr  -  /  Kit)  n  RUj)  exp.  -/(  R(u)duJdP 


cN(t) 


NUj) 


H  ’ t  N l x  )  exp(-«  l  (x  -T 
1-1  j-1  1  J 


DjdPfSfXj).  •••  ,N(  t)  ,T) 


N(x  ) 

•  n  n  i 

E!«n  (  II  Six  ))N(t)  exp(-  <  }  £  lx  — T  ))  . 

1-1  1  i-1  1-1  1  1 


(2.11 • 


Tills  expression  is  analogous  to  (2.8).  Hence,  relabeling  the  x^'s  and  t 
as  v  '»  as  in  (2.9),  it  is  easily  shown  that  (2.11)  is  equal  to 

n-*-! 

En<>i>- •"  •N<!'n*i>ii;iiS(>'i>'kr  -•  •’i(W'Vi:<  'j-Vi’ 

n-i 

*  exp  - <  i  l  (y  -T  ) : ]) 

j-l.J^i  n-1 


n+l 

•  t  exp ( - vy 


«r  kn+l  ^2  n+l  k 


j.VV'j-’ 


y  )  l  I  ...  I  UUv" 

kn.l-°V°  kl'°  J’1  ’ 

k.-k  i-1  0  n+l 

j>  *  exp.  t  l  ( (n-J+1  )y  -  J  y^ 

J-1  m-j 


n+1  n+1 

*  <ki-k4  l*  4  f  l  Wn- J+2)y  -  l  y  )<k  -k’  ) 1 

1-i+l  J  m-1  ' 


/  l  <k  -k  >!  . 
J-1  J  J 


x  j  <  t  s  *1  ,  i  -  1,  •  •  •  ,n  ♦  1  , 


(2.12) 


where  y^+  ^  t  ( y  t  ^  ,yi  +  1  as  in  (2.6).  k';  •  k^_j  lor  J  s  1  ,  and 

k'  -  k^  for  J  >  1  . 

Therefore,  the  haves  estimator  K(t)  of  r(t)  is  given  by  (2.12) 
divided  by  (2.7).  This  estimator  is  a  smooth  curve,  rather  than  a  step 
function  as  the  estimator  of  .Marshall  and  Proscnan  (1965)  yields,  and  is  an 
Increasing  function  over  the  range  of  the  data  so  that  we  have  a  ''closure’’ 
property  of  I  FR  over  this  range. 

We  remark  that  the  constant  Jump  size  t  mav  be  replaced  by  a  random 
jump  size  W  ,  sav,  with  distribution  function  (i  .  This  extension  compli¬ 
cates  the  expressions  and  makes  them  intractable.  The  constant  Jump  prior 
process  seems  reasonable  and  sufficient  for  Bavesian  analysis  as  was  dis¬ 
cussed  in  Section  1. 


3.  AS  EXAMPLE 

As  an  example,  the  formulas  for  the  Bavesian  nonpa raise  trie  estimators 
derived  in  Section  2  were  programmed  for  computation  on  an  electronic  com¬ 
puter.  A  random  sample  of  size  n  *  5  was  selected  fr  -m  the  100  observa¬ 
tions  of  Dubey  (1967)  generated  from  a  Weibull  distribution  with  scale 
parameter  one  and  shape  parameter  1.2,  so  that  the  failure  rate  function  is 
increasing.  The  sample  selected  given  in  order  is  0.135873,  0.666654, 
0.948871,  1.341265,  1.521437.  Figures  1  and  2  show  the  resulting  estimates 
of  r(t)  and  F(t)  ,  respectively.  Several  values  of  the  prior  process 
parameters  v  and  <  were  used  to  illustrate  the  effects  of  the  choices. 
Also,  computed  and  plotted  in  the  figures  were  the  maximum  likelihood  esti¬ 
mates  (Marshall  and  Prcschan  (1/65))  and  the  true  r(t)  and  F(t)  for 
comparison.  The  mle  of  r(t)  is  ♦  “P  at  1.521437  and  beyond. 


I m.it *•««  of  F-illtir 
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The  computations  in  Figures  1  and  2  were  performed  with  13  terms  in  the 
outside  infinite  summations  in  expressions  (2.7),  (2.10),  and  (2.12).  The 

results  were  t  e  sane  to  at  least  two  decimal  places  when  using  ten  or 

(  0 

nineteen  terms.  Also,  an  approximation  to  thi  x  in  (2.b)  (or  in 

0  xi  xi-l 

the  numerator  formulas)  given  by  the  midpoint  x  m  - - -  yielded  the 


sane  results  to  at  least  two  decimal  places. 


THf  CAST  OF  AKB  1  IRAK  1  l.V  R 1 CHT - CLNSORED  DATA 


In  many  situations  in  lift-  testing  and  survival  analysis  the  items  mav 

be  entered  into  and  removed  from  the  studv  at  arbitrary  times.  Susarla  and 

Van  Ryxin  (197b)  and  Ferguson  and  Fhadia  (1979)  havt  obtained  Bayesian  non- 

paranetrlc  cstirutors  of  survival  probability  utilizing  the  Dirichlet 

processes  and  processes  neutral  to  the  right  as  priors,  respectively. 

Let  X'.  •••  ,xs  denote  a  random  sample  of  true  survival  times  of  n 
1  n 

items  whose  lifetime  distribution  is  F  .  Let  l',,  •••  .1'  be  a  set  of 

J  n 

constants  or  independent  random  variables  which  are  also  independent  of 

X*  ,  • • •  ,  X*  .  Define 
1  n 


X.  •  min {X*,L  j )  . 


P 

1 

0 


if  X, 


if  X*  >  V 


i  •  1 ,  • •  •  ,n  . 


The  pairs  (X^.'^i  •  1  "  1,  *••  ,n  ,  are  the  observat ions  which  are  available 
to  the  experimenter,  in  which  X^  is  an  observed  lifetime  if  1^*1  and  a 
censored  lifetime  if  6^  -  0  .  That  is,  it  Is  known  which  observations 
represent  ''failures''  and  which  ones  represent  "losses". 
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The  likelihood  function  can  be  written  as 


L  •  II  [f(x  >]  i[F(x  ) 
1-1 


-  H  lR(x  )]  1  Fix  ) 

t  _  1  *  * 


i Lagakos 


1-6 


i 


(1979)) 


where  R  ( t  J  is  the  failure  rate  function  as  before.  Hence,  for  the  right- 
censored  sample  (X  , ' j )  ,  i  •  1,  •••  ,n  ,  the  Bayesian  nonparanet r lc  esti¬ 
mator  of  F(t)  is  obtained  as  in  Section  2  for  complete  samples.  Die 
estimator  is  given  h\ 


F  <  t )  •  1  F(t)|(xt,6  ),i-l,  •••  ,n 

_  n  xi 

Fit'  Kix  )  exp [-J  Riu)du  "dP 

' _ _  i-1 _ _ _ ° _ 

n  (  Xj 

!R(x  )  exp!-/  R(u)du]dP 
0  i-1  0 


(4.1  ' 


whe  Ti¬ 


the  notation  .  the  sartu  as  before. 

The  Denominator  of  (4.1)  i .»  given  by 

k  k, 

<*•  n  2  n 


cxp(-vx  )  i  l  •••  [  H  (*k 

n  k  -0  k  -0  k  -0  1-1 

n  n- 1  1 


with 


X,  S  •  •  •  S  X 
1  n 


1  •  1 


■  exp 


i(n-J'l)x  -  l  x  )(k  -k  >1/  T  (k  -k  )!  . 
)-l  '  i-J  5  1-1 


(4.2) 


For  the  numerator  of  (4.1),  we  relabel 

similar  to  that  in  Section  2.  When  », 

x_  •  0  and  x  ,  -  -  ,  let 
0  n*l 


i-1 


the  observations  In  a  manner 
t  t  x  ,  i-1,  •••  ,n  ♦  1  ,  with 


0  and  x 
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yj  ‘  V'3 


xj  •  j  -  i  -  1 
t  . 

xj-l 


J  -  1 

.  j  »  i  ♦  1 


and  Y 


.  J  «  i  -  1 

0  .  j  -  1 

Vi  •  3  * 


Then  the  numerator  is  written  analogous  to  (2.10)  at 

3  *  n+l 


kn-H  k2  n+1 


exp(-v>n+i )  l  l 


l  l  •••  l  t!  (  <  k  )  J  V  n+1/  n  <k  -k  )! 

k  k  -0  k  -0  1-1  5  1-1  ^ 

n+1-0  n  1  J  3 


n+1  krki-i 

"  Oj-Vj.,)  3  3 

j-1  J  J  * 


n+1  n+1 

exp  •  l  ( (n- J+l )y .  -  l  V  >  (k  -k  . )  . 

J-1  3  n>-  J  '  3  3 


xi-l  <  1  =*  xi  '  1  "  *n4  1 


(4.3) 


Therefore,  for  right-censored  samples  the  haves  nonparametr ic  estimator 
F(l>  of  F(t)  Is  given  bv  (-.3)  divided  bv  (-4.2).  The  nurerator  foi  the 
haves  estimate  of  r(t)  can  be  found  in  a  similar  manner  and  is  omitted  here. 

Remark  1.  We  remark  that  the  havesian  nonparametric  estimator  for 
F(t)  given  here  is  a  smooth  curve  (and  assumes  F  is  an  IFR  distribution). 
Susarla  and  Van  Ryzin’s  (1976)  havesian  nonparametric  estimator  for  this  case 
using  the  Dirichlet  priors  is  discontinuous  at  the  uncensored  data  points 
while  the  well-known  F.aplan  and  Meier's  (19S8)  estimator  is  a  step  function 
with  jumps  at  the  uncensored  observations. 

Remark  2.  It  should  also  be  noted  that  we  can  obtain  a  Bayesian  nonparametric 

estimator  of  the  densitv  function  f(t)  in  this  framework.  For  the  censored 

sample  1  ”  * . n,  w‘‘  calculate  f  n  ( t )  -  E'f(t)  |  (x^,5^),  *  "  ^ . n 

in  a  manner  similar  to  that  for  (4.1),  Since  f(t)  -  R(t)  F(t),  the  required 

numerator  for  f  (t)  is  the  same  as  (4.3)  with  the  stated  relabeling  of  x  '»  and 
n  i 


u 


t  to  v , ' K  dn<*  -V's  to  >,'s  except  that  ■  1  for  j  ■  i.  Again,  fn(t)  *b  a 
smooth  density  estimator  and,  of  course,  for  6^  “  1  for  all  i,  we  have  the 
case  of  no  censoring  as  considered  in  Section  2. 
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